Transition from pure-state to mixed-state entanglement by random scattering 
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We calculate the effect of polarization-dependent scattering by disorder on the degree of 
polarization-entanglement of two beams of radiation. Multi-mode detection converts an initially 
pure state into a mixed state with respect to the polarization degrees of freedom. The degree of 
entanglement decays exponentially with the number of detected modes if the scattering mixes the 
polarization directions and algebraically if it does not. 
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I. INTRODUCTION 



A pair of photons in the Bell state (|HV) + |VH)) /\/2 
can be transported over long distances with little degra- 
dation of the entanglement of their horizontal (H) and 
vertical (V) polarizations. Polarization-dependent scat- 
tering has little effect on the degree of entanglement, as 
long as it remains linear (hence describable by a scatter- 
ing matrix) and as long as the photons are detected in 
a single spatial mode only. This robustness of photon 
entanglement was demonstrated dramatically in a recent 
experiment and theory 0, on plasmon-assisted en- 
tanglement transfer. 



Polarization-dependent scattering may significantly 
degrade the entanglement in the case of multi-mode de- 
tection. Upon summation over N spatial modes the ini- 
tially pure state of the Bell pair is reduced to a mixed 
state with respect to the polarization degrees of freedom. 
This loss of purity diminishes the entanglement — even 
if the two polarization directions are not mixed by the 
scattering. 



The transition from pure-state to mixed-state entan- 
glement will in general depend on the detailed form of 
the scattering matrix. However, a universal regime is 
entered in the case of randomly located scattering cen- 
tra. This is the regime of applicability of random-matrix 
theory As we will show in this paper, the trans- 

mission of polarization-entangled radiation through dis- 
ordered media reduces the degree of entanglement in a 
way which, on average, depends only on the number N 
of detected modes. (The average refers to an ensemble of 
disordered media with different random positions of the 
scatterers.) The degree of entanglement (as quantified 
either by the concurrence [|| or by the violation of a Bell 
inequality 0, H|) decreases exponentially with N if the 
disorder randomly mixes the polarization directions. If 
the polarization is conserved, then the decrease is a power 
law (cx TV -1 if both photons are scattered and cx TV -1 / 2 
if only one photon is scattered). 
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FIG. 1: Schematic diagram of the transfer of polarization- 
entangled radiation through two disordered media. The de- 
gree of entanglement of the transmitted radiation is measured 
by two multi-mode photodetectors (Ni modes) in a coinci- 
dence circuit. The combination of polarization-dependent 
scattering and multi-mode detection causes a transition from 
a pure state to a mixed state in the polarization degrees of 
freedom, and a resulting decrease of the detected entangle- 
ment. 



II. FORMULATION OF THE PROBLEM 

We consider two beams of polarization-entangled pho- 
tons (Bell pairs) that are scattered by two separate dis- 
ordered media (see Fig.^l. Two photodetectors in a co- 
incidence circuit measure the degree of entanglement of 
the transmitted radiation through the violation of a Bell 
inequality. The scattered Bell pair is in the pure state 
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( U na V nir + U uo V ™t) 



(2.1) 



The indices n € {1,2,... Mi}, m G {1,2,... M 2 } la- 
bel the transverse spatial modes and the indices a, r G 
{+,— } label the horizontal and vertical polarizations. 
The first pair of indices n, a refers to the first photon 
and the second pair of indices m, t refers to the second 
photon. The scattering amplitudes u+ a , u~ a of the first 
medium form two orthonormal vectors of length M\, and 
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similarly for the second medium. 

A subset of Ni out of the Mi modes are detected 
in the first detector. We relabel the modes so that 
n = 1,2, ... Ni are the detected modes. This subset is 
contained in the four vectors w+ + = u^ + , u+~ = u^_, 
u~ + = u~ + , u~~ = u~_ of length Ni each. We write 
these vectors in bold face, u±±, omitting the mode in- 
dex. Similarly, the second detector detects N2 modes, 
contained in vectors v±±. A single or double dot be- 



tween two pairs of vectors denotes a single or double 
contraction over the mode indices: a • b = ^2 n a n b n , 
ab : cd = J2 n ,m a n b m c m d n - 

The pure state has density matrix ^ na mr ^* nl(J , m i T i- 
By tracing over the detected modes the pure state is re- 
duced to a mixed state with respect to the polarization 
degrees of freedom. The reduced density matrix is 4 x 4, 
with elements 



Par.a'T' = (u +(T V_ r + U_ ff V +r ) I (vl T ,U^ CT , + V^U*^), (2.2) 

Z = J2 (u +CT v_ r + u CT v +r ) : {v*_ T u* +a + v* +T u*_ a ) . (2.3) 

<T,T 

The complex numbers that enter into the density matrix are conveniently grouped into a pair of Hermitian positive 
definite matrices a and b, with elements cicrT,(T'T' — u aT ' u a'r'> ^ot.o't' = v CTT • v*, r ,. One has 

ZP(JT,<J'T' = 0- + <J,+<j'b-T,-T' + <l-a,-Cr'b+T,+T' + (l-B,+B'b+T,— t' + a +<7, -IT' b— T,+t' ■ (2-4) 

I 



The degree of entanglement of the mixed state with 
4x4 density matrix p is quantified by the concurrence 
C, given by 

C = max|o,V^i- V^2-V^3- V^}- ( 2 - 5 ) 

The Xi 's are the eigenvalues of the matrix product 

p- (crj, ® a y ) ■ p* ■ (a y ® a y ), 

in the order Ai > A2 > A3 > A4, with a y a Pauli matrix. 
The concurrence ranges from (no entanglement) to 1 
(maximal entanglement) . 

In a typical experiment , the photodetectors can not 
measure C directly, but instead infer the degree of en- 
tanglement through the maximal violation of the Bell- 
CHSH (Clauser-Horne-Shimony-Holt) inequality 0, @. 
The maximal value £ of the Bell-CHSH parameter for 
an arbitrary mixed state was analyzed in Refs. 0, . 
For a pure state with concurrence C one has simply 
For a mixed state there is no one-to- 
one relation between £ and C. Depending on the den- 
sity matrix, £ can take on values between 2Cy/2 and 
2Vl + C 2 , so £ > 2 implies C > but not the other 
way around. The general formula 

£ = 2 Vmi + u 2 (2.6) 

for the dependence of £ on p involves the two largest 
eigenvalues u\,U2 of the real symmetric 3x3 matrix R T R 
constructed from Rki = Tr p Uk <S> 07 . Here a± , 02 , 03 refer 
to the three Pauli matrices a x ,a y , a z , respectively. 

We will calculate both the true concurrence C and the 
pseudo-concurrence 

C = Vmax(0, £ 2 /A- 1) < C (2.7) 



inferred from the Bell inequality violation. 

As a special case we will also consider what happens if 
only one of the two beams is scattered. The other beam 
reaches the photodetector without changing its mode or 
polarization, so we may set v^ lrT = S m: iS a .±. This implies 
baT.a'T' — ^a,T^a / ,r' ; hence 

where we have defined f = — r. The normalization is now 
given simply by Z = J2 a , T a <rr,aT- 



III. RANDOM-MATRIX THEORY 

For a statistical description we use results from the 
random-matrix theory (RMT) of scattering by disordered 
media 0,0. According to that theory, the real and imag- 
inary parts of the complex scattering amplitudes are 
statistically distributed as independent random variables 
with the same Gaussian distribution of zero mean. The 
variance of the Gaussian drops out of the density matrix; 
we fix it at 1. The assumption of independent variables 
ignores the orthonormality constraint of the vectors u, 
which is justified if TVi <§; Mi. Similarly, for N2 <C M2 
the real and imaginary parts of vfj have independent 
Gaussian distributions with zero mean and a variance 
which we may set at 1. 

The reduced density matrix of the mixed state depends 
on the two independent random matrices a and b, ac- 
cording to Eq. H2.4f) . The matrix elements are not inde- 
pendent. We calculate the joint probability distribution 
of the matrix elements, using the following result from 
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RMT ^| : Let W be a rectangular matrix of dimension 
px (k+p), filled with complex numbers with distribution 

P({W nm }) oc exp (-cTrWW 1 ") , c> 0. (3.1) 

Then the square matrix H = WW^ (of dimension p x p) 
has the Laguerre distribution 

P{{H nm }) cx (Dctff) fc cxp(-cTr£0. (3.2) 

Note that H is Hermitian and positive definite, so its 
eigenvalues h n (n = 1, 2, . . .p) are real positive numbers. 
Their joint distribution is that of the Laguerre unitary 
ensemble, 

P({hn}) cx J] h k n e~ ch " Hihi - h 3 f. (3.3) 

n i<j 

The factor (hi — hj) 2 is the Jacobian of the transforma- 
tion from complex matrix elements to real eigenvalues. 
The eigenvectors of H form a unitary matrix U which is 
uniformly distributed in the unitary group. 

To apply this to the matrix a we set c = 1/2, p = 4, 
k = Ni — 4. We first assume that N% > 4, to ensure that 
k > 0. Then 

P({a aT ^ T ,}) cx (Deta) Ari ~ 4 exp(-iTra) , (3.4) 

P({a n }) oc J] a^ 4 e~ Q " /2 J]^ " a i)^ (3-5) 

where oi, 0,2, as, a 4 are the real positive eigenvalues of a. 
The 4x4 matrix U of eigenvectors of a is uniformly 
distributed in the unitary group. If Ni = 1,2,3 we set 
c = 1/2, p = Ni, k = A-N x . The matrix a has 4 - N x 
eigenvalues equal to 0. The N\ non-zero eigenvalues have 
distribution 

P({a n }) « J] atr Nl e~ a "/ 2 J]^ " ■ (3-6) 

n i<j 

The distribution of the matrix elements b aT ,a'T' and of 
the eigenvalues b n is obtained upon replacement of Ni 
by N 2 in Eqs. $TQ . J23J), and (|33|) . 

IV. ASYMPTOTIC ANALYSIS 

We wish to average the concurrence (|2.5f) and pseudo- 
concurrence (|2.7|l with the RMT distribution of Sec. ITTT1 
The result depends only on the number of detected modes 
N\, N2 in the two photodetectors. Microscopic details of 
the scattering media become irrelevant once we assume 
random scattering. The averages (C), (C) can be cal- 
culated by numerical integration |l3j |. Before present- 
ing these results, we analyze the asymptotic behavior 
for JV; > 1 analytically. We assume for simplicity that 
Ni=N 2 = N. 

It is convenient to scale the eigenvalues as 

a n = 27V(l + a„), b n = 2N{\ + f3 n ). (4.1) 



The distribution of the a„'s and /3„'s takes the same form 

P({a n }) cx exp ^-AT^[a„-ln(l + a n )] + 0(l)J , 

(4.2) 

where 0(1) denotes TV-independent terms. The bulk of 
the distribution l|4.2() lies in the region J2 n a n ^ 1/TV <C 
1, localized at the origin. Outside of this region the dis- 
tribution decays exponentially cx exp[— N f({a n })], with 

4 

/({««}) = ^K-ln(l + «„)]. (4.3) 

n=l 

The concurrence C and pseudo-concurrence C depend 
on the rescaled eigenvalues a n ,(3 n and also on the pair 
of 4 x 4 unitary matrices U, V of eigenvectors of a and b. 
Both quantities are independent of N, because the scale 
factor N in Eq. Ij4.1|l drops out of the density matrix 
(|2.4(l upon normalization. 

The two quantities C and C are identically zero when 
the a n 's and /3„'s are all <C 1 in absolute value. For 
a nonzero value one has to go deep into the tail of the 
eigenvalue distribution. The average of C is dominated 
by the "optimal fluctuation" a° pt , (3° pt , U opt , V opt of 
eigenvalues and eigenvectors, which minimizes f({cv n }) + 
f{{Pn}) in the region C > 0. The decay 

(C) ^ exp {-N[f{{aT }) + /({/?° pt })]) = e~ AN (4.4) 

of the average concurrence is exponential in TV, with a 
coefficient A of order unity determined by the optimal 
fluctuation. The average (C) ~ e~ BN also decays ex- 
ponentially with N, but with a different coefficient B in 
the exponent. The numbers A and B can be calculated 
analytically for the case that only one of the two beams 
is scattered. 

Scattering of a single beam corresponds to a density 
matrix p which is directly given by the matrix a, cf. Eq. 
1)2.8(1 . To find A, we therefore need to minimize f{{a n }) 
over the eigenvalues and eigenvectors of a with the con- 
straint C > 0, 

A = rain rr {/({«„}) | C (p({a n },U)) > 0} . (4.5) 

The minimum can be found with the help of the following 
result The concurrence C (p) of the two-qubit density 
matrix p, with fixed eigenvalues Ai > A2 > A3 > A4 
but arbitrary eigenvectors, is maximized upon unitary 
transformation by 

maxC(f7pfi t ) = max jo, A x - A 3 - 2 VA 2 A 4 | . (4.6) 

(The matrix fl varies over all 4 x 4 unitary matrices.) 
With this knowledge, Eq. I|4.5|l reduces to 

A = min {f({a n })\ ai - a 3 - 2 y/(l + a 2 )(l + a 4 ) > 0}, 

{a„} 

(4.7) 



4 



where we have ordered a± > a 2 > > This yields 
for the optimal fluctuation a° pt = 1, a° pt = a^ pt = 
^opt _ _j^3 anc [ 

,4 = 31n3-41n2 = 0.523. (4.8) 
The asymptotic decay (C) oc e~ AN is in good agreement 
I 



hence 

S = ln(ll + 5V5)-ln2 = 2.406. (4.12) 

The decay (C) oc e~ BN is again in good agreement with 
the numerical results for finite N (Fig. |2J. 

If both beams are scattered, a calculation of the opti- 
mal fluctuation is more complicated because the eigen- 
values {a n }, {/?„} and the eigenvectors U, V get mixed 
in the density matrix l|2.4|l . The numerics of Fig. [21 gives 
(C) cx e - 3 3JV for the asymptotic decay of the concur- 
rence. The averaged pseudo-concurrence for two-beam 
scattering could not be determined accurately enough to 
extract a reliable value for the decay constant. 



V. COMPARISON WITH THE CASE OF 
POLARIZATION-CONSERVING SCATTERING 

If the scatterers are translationally invariant in one di- 
rection, then the two polarizations are not mixed by the 
scattering. Such scatterers have been realized as paral- 
lel glass fibers [l5|. One polarization corresponds to the 
electric field parallel to the scatterers (TE polarization), 
the other to parallel magnetic field (TM polarization). 
The boundary condition differs for the two polarizations 
(Dirichlet for TE and Neumann for TM), so the scatter- 
ing amplitudes u ++ , v ++ , u , v that conserve the 

polarization can still be considered to be independent 
random numbers. The amplitudes that couple different 
polarizations vanish: , v_| , u |_, v |_ are all zero. 

The reduced density matrix l|2.4|) simplifies to 



with a numerical calculation for finite N, see Fig. |5J 



The asymptotic decay of the average pseudo- 
concurrence (C) for a single scattered beam can be found 
in a similar way, using the result [To| 



(4.9) 



(4.10) 



(4.11) 
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FIG. 2: Average concurrence (C) (squares) and pseudo- 
concurrence (C'} (triangles) as a function of the number N 
of detected modes. Closed symbols are for the case that only 
one of the two beams is scattered and open symbols for the 
case that both beams are scattered. The decay of (C') in 
the latter case could not be determined accurately enough 
and is therefore omitted from the plot. The solid lines are 
the analytically obtained exponential decays, with constants 
A = 3 In 3 - 4 In 2 and B = ln(ll + 5\/5) - In 2, cf. Eqs. 
and (£13 • 



with r = — t, t' = — t'. We will abbreviate A aT = 
o-oo,tt, B aT = b aa _ TT . The concurrence C and pseudo- 



maxC^pQt) = Vmax {0, 2[k x - A 4 ) 2 + 2(A 2 - A 3 ) 2 - (Ai + A 2 + A 3 + A 4 ) 2 }. 

To obtain the optimal fluctuation we have to solve 

B = min {/({a„})|2(ai - a 4 ) 2 + 2(a 2 - a 3 ) 2 - (4 + a Y + a 2 + a 3 + a A f > 0}, 

which gives 

a° pt = i(-l + 2V2 + V5), aT = aT = \{l-^>), a%* = ±(-1- 2y/2 + y/E), 

I 
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concurrence C are calculated from Eqs. (|2.5|l and l|2.7|l . 
with the result 



C = C 



2\A+-\\B+- 



(5.2) 



++ 



It is again our objective to calculate (C) for the case 
N% = N2 — N. The distribution of the matrices A and 
B follows by substituting N\ — 4 — > N — 2 in Eq. 



P{{A aT }) cx (DetA) JV - 2 eocp(-|'&A) 



(5.3) 



The average over this distribution was done numerically, 
see Fig. |3J For large N we may perform the following 
asymptotic analysis. 

We scale the matrices A and B as 



A = 2N{±+A), B = 2N(t + B). 



(5.4) 



In the limit N — > 00 the Hcrmitian matrices .A and 23 
have the Gaussian distribution 



^({•Arr}) 



(5.5) 



(The same distribution holds for B.) In contrast to the 
analysis in Sec. II VI the concurrence does not vanish in the 
bulk of the distribution. The average of Eq. 15.211 with 
distribution H5.5JI yields the algebraic decay 



<0 = ~, ">1, 



(5.6) 



in good agreement with the numerical calculation for fi- 
nite N (Fig.EJ. 

A completely analytical calculation for any N can be 
done in the case that only one of the beams is scattered. 
In that case B aT — 1 and the concurrence reduces to 



C 



2|A+_| 



.4 



++ 



A— 



(5.7) 



Averaging Eq. (|5.7|l over the Laguerre distribution (|5.3|) 
gives 



(C) 



/ttT(N + 1/2) 



(5.8) 



2 r(AT+l) ' 
For large N, the average concurrence (|5.8|l falls off as 

(C) = ^-±=, N»l. (5.9) 

This case is also included in Fig. 

VI. CONCLUSION 

In summary, we have applied the method of random- 
matrix theory (RMT) to the problem of entanglement 
transfer through a random medium. RMT has been used 
before to study production of entanglement |l6l Il7l ITsl 
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FIG. 3: Average concurrence (C) as a function of the number 
N of detected modes, for the case of polarization-conserving 
scattering of both beams (open squares) and one beam (closed 
squares). The data points are the result of a numerical av- 
erage. The dashed line is the asymptotic result 15.61 and 
the dotted line is the analytical result 15.81 . The pseudo- 
concurrence C' is identical to C for polarization-conserving 
scattering. 



El HI HI El E3. Here we have studied the loss of 
entanglement in the transition from a pure state to a 
mixed state. 

A common feature of all these theories is that the re- 
sults are universal, independent of microscopic details. 
In our problem the decay of the degree of entanglement 
depends on the number of detected modes but not on 
microscopic parameters such as the scattering mean free 
path. 

The origin of this universality is the central limit the- 
orem: The complex scattering amplitude from one mode 
in the source to one mode in the detector is the sum 
over a large number of complex partial amplitudes, cor- 
responding to different sequences of multiple scattering. 
The probability distribution of the sum becomes a Gaus- 
sian with zero mean (because the random phases of the 
partial amplitudes average out to zero). The variance 
of the Gaussian will depend on the mean free path, but 
it drops out upon normalization of the reduced density 
matrix. The applicability of the central limit theorem 
only requires that the separation of source and detector 
is large compared to the scattering mean free path, to 
ensure a large number of terms in the sum over partial 
amplitudes. 

The degree of entanglement (as quantified by the con- 
currence or violation of the Bell inequality) then de- 
pends only on the number N of detected modes. We 
have identified two qualitatively different types of de- 
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cay. The decay is exponential oc e~ if the scattering 
mixes spatial modes as well as polarization directions. 
The coefficient c depends on which measure of entangle- 
ment one uses (concurrence or violation of Bell inequal- 
ity) and it also depends on whether both photons in the 
Bell pair are scattered or only one of them is. For this 
latter case of single-beam scattering, the coefficients c 
are 3 In 3 — 4 In 2 (concurrence) and ln(ll + 5\/5) — In 2 
(pseudo-concurrence). The decay is algebraic oc N~ p if 
the scattering preserves the polarization. The power p 
is 1 if both photons are scattered and 1/2 if only one of 
them is. Polarization-conserving scattering is special; it 
would require translational invariance of the scatterers in 
one direction. The generic decay is therefore exponential. 
Finally, we remark that the results presented here ap- 



ply not only to scattering by disorder, but also to scat- 
tering by a cavity with a chaotic phase space. An experi- 
mental search for entanglement loss by chaotic scattering 
has been reported by Woerdman et al. |24j . 
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